function res = GroupRbSpinPolarization(obj, aomVoltageList, pumpPowerList, pumpFrequency, probePower, probeFrequency, cellTemperature, varargin)
    p=inputParser;
    p.addRequired('aomVoltageList', @isnumeric);
    p.addRequired('pumpPowerList', @isnumeric);
    p.addRequired('pumpFrequency', @isnumeric);
    p.addRequired('probePower', @isnumeric);
    p.addRequired('probeFrequency', @isnumeric);
    p.addRequired('cellTemperature', @isnumeric);
    
    p.addParameter('resonant_frequency', 100e3, @(x) x>0);
    p.addParameter('driving_amplitude', 0.1, @(x) x>0);
    p.addParameter('guessField', 6.0, @(x) x>0);
    p.addParameter('fittingWidth', 1.0, @(x) x>0);
    p.addParameter('bz_min', 0.0, @(x) x>=0);
    p.addParameter('bz_max', 7.0, @(x) x>=0);
    p.addParameter('npoint', 100, @(x) x>0);
    p.addParameter('average', 100, @(x) x>0);
    p.addParameter('order', 4, @(x) ismember(x, 1:8));
    p.addParameter('bandwidth', 100, @(x) x>0);
    p.addParameter('description', 'normal', @ischar);
    p.parse(aomVoltageList, pumpPowerList, pumpFrequency, probePower, probeFrequency, cellTemperature, varargin{:});

    %%
    nV = length(p.Results.aomVoltageList);
    if length(p.Results.pumpPowerList) ~= nV
        error('Pump Power must be the same length of Vaom.');
    end
    
    res.vAOM = p.Results.aomVoltageList;
    res.pumpPower = p.Results.pumpPowerList;
    res.data = cell(1, nV);
    res.width = zeros(1, nV);
    res.hight = zeros(1, nV);
    res.bz0 = zeros(1, nV);
    %%
    getReporter;
    rpt.NewGroup(sprintf('Rb Spin Polarization vs Pump Power, T=%d degC, ProbePower = %5.3f', cellTemperature, probePower));
    rpt.MainText(['In this group, we measure RF resonance peaks under pumping with different power Vaom = ', sprintf('%fV, ', aomVoltageList)]);
    rpt.TextEnter();
    
    dataCell = {'pump freq (THz)',                         'probe power (mW)',                     'probe freq (THz)',                         'cell temp (degC)',                          'Rb freq (Hz)',                                 'drv amp (V)'; ...
                sprintf('%8.5f', p.Results.pumpFrequency), sprintf('%3.2f', p.Results.probePower), sprintf('%8.5f', p.Results.probeFrequency), sprintf('%3.1f', p.Results.cellTemperature), sprintf('%3.1e', p.Results.resonant_frequency), sprintf('%5.3f', p.Results.driving_amplitude)};
    rpt.AppendTable(dataCell, 'Unchanged Parameters.');
    res.condition = dataCell;
    
    for k = 1:nV
        vAOM = p.Results.aomVoltageList(k);
        if vAOM == 0
            obj.aom.disable(); pause(1.0);
            res_data = obj.GetMainField('description', sprintf('AOM Output disabled'), 'isReport', true);
        else
            obj.aom.setDC(vAOM); pause(1.0);
            res_data = obj.GetMainField('description', sprintf('vAOM=%2.1f V', vAOM), 'isReport', true);
        end
        res.data{k}  = res_data;
        res.width(k) = res_data.width;
        res.hight(k) = res_data.hight;
        res.bz0(k) = res_data.bz0;
    end
    
    fig = obj.PlotRbSpinPolarization(res, p.Results.resonant_frequency/p.Results.guessField);
    rpt.NewRecord('Summary');
    rpt.AppendFigure(fig, 'Summary');
    
    getSession;
    sess.addData('RbPolarizationGroup', res, p.Results.description);
end

